Bayesian Linear Analysis

Changes since Meeting - Tested covariance of station 51N; included in GoA not northern CC - PDO is unstandardized - Upwelling is standardized across entire time period (1967 - 2022), within region and within season (for the seasonal models - annaual model is not) - Tested model using seasonal averages instead of monthly means across season - Added a lagged year variable so the months of Nov and Dec are inculded with the correct Jan - Mar values for winter

Mapping Bakun Index locations with the regions we have identified

map

Combining Bakun Data

1. Import the Bakun upwelling data from each individual file associated with the 13 lat/long locations. Data versions are: 1˚ 6-hourly version

2. Extract yyyy-mm-ddThh:mm:ss from the ‘time’ column and add region indicator based on station ID and add a time period based on year

## `summarise()` has grouped output by 'station_id', 'Year', 'YearDay'. You can
## override using the `.groups` argument.
## Joining with `by = join_by(period)`

## Warning: Removed 1365 rows containing missing values (`geom_line()`).

## # A tibble: 1 × 7
##   station_id  Year region upwelling_index_cum YearDay period period2    
##   <chr>      <dbl> <chr>                <dbl>   <dbl> <chr>  <chr>      
## 1 57N         1968 GoA                -45081.     366 1      1967 - 1988

## `geom_smooth()` using method = 'gam' and formula = 'y ~ s(x, bs = "cs")'

## Warning: Removed 17745 rows containing non-finite values (`stat_smooth()`).
## Removed 1365 rows containing missing values (`geom_line()`).

2. Test covariation and correlation between station 51N and 48N/54N. Station 51N covaries more with both stations 54N and 57N (GoA) compared to 48N (Northern CC) and I elected to include it in the GoA

3. Summarize across days and hours to get monthly means for region and year. If you revert to seasonal summaries return to here to summarise across dates without taking means of means. Adding standardized upwelling index - upwelling is standardized across the entirety of the time period, but within region and space. Note because of the standardization positive values are above average for a given region and season for a year (cannot be interpreted as downwelling versus upwelling)

4. Import PDO and NPGO data. Label months by a season factor. NOT standardized in this version - Note PDO CAN be interpreted as positive or negative BUT this is the data from the website so it has been detrended.

24. Merging the standardized datasets into a single dataframe and assigning levels

Run the STAN model outside of markdown file using “STANrun.R” which exports the posteriors so that the model does not run in markdown when you knit.

STAN Model

Bayesian linear model with era specific intercept and slope and region specific slope,

\[Upwelling = \alpha_{e,r} + \beta_{e,r}*PDO +\sigma\] where e is factor “era” corresponding to eras 1967 - 1988, 1989 - 2013, and 2014 - 2022 and r corresponds to regions Southern California Current (south of Mendocino), Northern California Current (Mendocino to Vancouver Island), and Gulf of Alaska (north of Vancover Island).

Parameter	Description	Prior
\(\alpha_{e,r}\)	Era and region specific anomaly	\(~ Normal(0,10)\)
\(\beta_{e,r}\)	Upwelling-PDO relationship by era and region	\(~ Normal(0,10)\)
\(\sigma\)	Combined observation and process error	\(~ Normal(0,10)[0,]\)

#PDO - Upwelling relationship analyses

Monthly Mean Model

Linear relationship between monthly PDO Index and monthly upwelling (Bakun 1˚ 6-hourly) for each era and each region estimated from a Bayesian regression model.

Era specific intercepts across the entire region (Southern California Current to Gulf of Alaska) estimated from a Bayesian regression model

Winter Model (Nov - March)

Using monthly mean data through winter months (Nov - Mar) where data is standardized for the time period and season

Linear relationship between monthly PDO Index and monthly upwelling (Bakun 1˚ 6-hourly) for each era and each region estimated from a Bayesian regression model using Winter Monthly Means (Nov-March).

Era specific intercepts across the entire region (Southern California Current to Gulf of Alaska) estimated from a Bayesian regression model for winter months (Nov - March)

Spring Model (April - June)

Using monthly mean data through winter months (April - June) where data is standardized for the time period and season

Linear relationship between monthly PDO Index and monthly upwelling (Bakun 1˚ 6-hourly) for each era and each region estimated from a Bayesian regression model using spring Monthly Means (April - June).

Era specific intercepts across the entire region (Southern California Current to Gulf of Alaska) estimated from a Bayesian regression model for spring monthly means (April - June)

Summer Model (July - August)

Using monthly mean data through winter months (July - August) where data is standardized for the time period and season

Linear relationship between monthly PDO Index and monthly upwelling (Bakun 1˚ 6-hourly) for each era and each region estimated from a Bayesian regression model using Summer Monthly Means (July - August).

Era specific intercepts across the entire region (Southern California Current to Gulf of Alaska) estimated from a Bayesian regression model for summer months (July - August)

Seasonal Mean Models

These results utilize seasonal averagers rather than monthly means within a season. May be better for identifying trends that are not spuriously seasonal (i.e. changes in a given month) but tradeoff is less data and may miss some interesting dynamics as a result

Annual Model

Plots of the annual average model.

## `geom_smooth()` using formula = 'y ~ x'

Winter Model

Fitting the winter seasonal model

Plots of the winter seasonal model.

## `geom_smooth()` using formula = 'y ~ x'

Spring Model

Fitting the spring seasonal model

Spring model plots

## `geom_smooth()` using formula = 'y ~ x'

## `geom_smooth()` using formula = 'y ~ x'

## Warning: The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?

Summer Model

Fitting the summer seasonal model

Summer seasonal plots

## `geom_smooth()` using formula = 'y ~ x'

#NPGO

Seasonal Mean Models

These results utilize seasonal averagers rather than monthly means within a season. May be better for identifying trends that are not spuriously seasonal (i.e. changes in a given month) but tradeoff is less data and may miss some interesting dynamics as a result

Annual Model

Plots of the annual average model.

## `geom_smooth()` using formula = 'y ~ x'

Winter Model

Fitting the winter seasonal model

Plots of the winter seasonal model.

## `geom_smooth()` using formula = 'y ~ x'

Spring Model

Fitting the spring seasonal model

Spring model plots

## `geom_smooth()` using formula = 'y ~ x'

## `geom_smooth()` using formula = 'y ~ x'

## Warning: The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?
## The following aesthetics were dropped during statistical transformation: label
## ℹ This can happen when ggplot fails to infer the correct grouping structure in
##   the data.
## ℹ Did you forget to specify a `group` aesthetic or to convert a numerical
##   variable into a factor?

Summer Model

Fitting the summer seasonal model

Summer seasonal plots

Questions / Thoughts

To Do:

1. Examine Upwelling relationships with SLP, SSH, wind stress, and SST Examine the relationship between SST and the same set of atmospheric variables (PDO, SLP, SSH, wind stress). Would we really expect different relationships with SLP given its reanalyses are what make up Bakun indices?

2. we could also compare these results to results with CUTI and BEUTI that use just the last two periods (CUTI and BEUTI account for more nuanced changes in cross and along shore wind; only data after 1988)

3. plankton case study I think that this paper would be more compelling if it had some degree of biology. Maybe we could bring in some of the Zoop CalCofi data as a case study within the larger regional analysis and link these changing relationships to an ecosystem response - it would also be a nice thing to build our way up the food web with the salmon/groundfish part of the project and still take a foodweb wide view that was initially in the proposal. Are there are plankton datasets in the Southern California Current and/or GoA we could use?

4. How should we think/talk about intercepts. To me slopes are more compelling…

Interpretations

NPGO

PDO

Winter across all regions, the 1967 - 1988 period shows the greatest difference in intercept with the 1989 - 2013 period where the most recent heat wave period actually falls between the two however slopes are much more consistent with substantial posterior overlap for all three periods in all regions.

Spring strong change in the slope of NCC during spring but less in other regions where it overlaps with the other two eras - this is consistent with both approaches to characterizing season.

The seasonal average model shows more of a difference for the SCC for the most recent heatwave period than the monthly mean model.

Summer Summer is pretty consistent in slop and intercept with the exception of the SCC intercept is different for the early and late periods compared to the middle.